{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Julia中的线性代数\n",
    "Based on work by Andreas Noack Jensen\n",
    "\n",
    "## 基础线性代数操作"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "首先定义一个随机矩阵"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = rand(1:4,3,3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "定义一个元素全为1的向量"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = fill(1.0, (3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "注意$A$的类型为Array{Int64,2}，而$x$的类型为Array{Float64,1}。Julia定义向量Vector{Type}的别名为Array{Type,1}，矩阵Matrix{Type}的别名为Array{Type,2}。\n",
    "\n",
    "许多基础操作和其他语言一样\n",
    "#### 乘法"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b = A*x"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 转置\n",
    "就像在其他语言中`A'`表示对`A`进行共轭转置"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们可以通过transpose获得转置矩阵"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "transpose(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 转置的乘法\n",
    "Julia中某些情况下可以省略`*`号"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A'A"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 解线性方程组\n",
    "用方阵$A$表示的线性方程组$Ax=b$用左除运算符（函数）`\\`求解"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A\\b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## 特殊的矩阵结构\n",
    "\n",
    "矩阵结构在线性代数中非常重要。接触一下大一些的线型系统就可以看到矩阵结构有*多*重要了。用线性代数标准包LinearAlgebra可以获得结构化的矩阵structured matrices："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "using LinearAlgebra"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 1000\n",
    "A = randn(n,n);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Julia可以推断特殊矩阵结构，比如判断对称矩阵"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Asym = A + A'\n",
    "issymmetric(Asym)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "但有时候浮点错误会比较麻烦"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Asym_noisy = copy(Asym)\n",
    "Asym_noisy[1,2] += 5eps()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "issymmetric(Asym_noisy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "幸运的是我们可以通过如`Diagonal`，`Triangular`，`Symmetric`，`Hermitian`，`Tridiagonal`和`SymTridiagonal`这样的函数来明确地创建矩阵"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Asym_explicit = Symmetric(Asym_noisy);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们来看看Julia计算`Asym`，`Asym_noisy`和`Asym_explicit`的特征值各要花多少时间"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "@time eigvals(Asym);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "@time eigvals(Asym_noisy);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "@time eigvals(Asym_explicit);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "本例中，使用`Symmetric()`处理`Asym_noisy` 后让计算效率提高了约5倍"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### A big problem\n",
    "Using the `Tridiagonal` and `SymTridiagonal` types to store tridiagonal matrices makes it possible to work with potentially very large tridiagonal problems. The following problem would not be possible to solve on a laptop if the matrix had to be stored as a (dense) `Matrix` type."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 1_000_000;\n",
    "A = SymTridiagonal(randn(n), randn(n-1));\n",
    "@time eigmax(A)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.0.5",
   "language": "julia",
   "name": "julia-1.0"
  },
  "language": "Julia",
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.0.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
